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Abstract 

We revisit the mean field model of globally and harmonically coupled 
parametric oscillators subject to periodic block pulses with initially random 
phases. The phase diagram of regions of collective parametric instability is 
presented, as is a detailed characterization of the motions underlying these 
instabilities. This presentation includes regimes not identified in earlier work 
[I. Bena and C. Van den Broeck, Europhys. Lett. 48, 498 (1999)]. In addi- 
tion to the familiar parametric instability of individual oscillators, two kinds 
of collective instabilities are identified. In one the mean amplitude diverges 
monotonically while in the other the divergence is oscillatory. The frequencies 
of collective oscillatory instabilities in general bear no simple relation to the 
eigenfrequencies of the individual oscillators nor to the frequency of the exter- 
nal modulation. Numerical simulations show that systems with only nearest 
neighbor coupling have collective instabilities similar to those of the mean 
field model. Many of the mean field results are already apparent in a simple 
dimer [M. Copelli and K. Lindenberg, to appear in Phys. Rev. E]. 
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I. INTRODUCTION 



Over the past two decades there has been an increasing interest in the nonequilibrium 
behavior of spatially extended systems modeled as ensembles of simple dynamical units cou- 
pled to each other. The collective evolution of such discrete coupled systems often exhibits 
qualitatively different behavior from that of the single units. An example of a system that 
has attracted a great deal of attention is a collection of a large number of coupled limit- 
cycle (phase) oscillators with randomly distributed natural frequencies PJ. This system 
has been invoked as a simple model for coupled chemical, biological or physical oscillators. 
A most spectacular collective phenomenon discovered with this model is a synchroniza- 
tion phase transition involving mutual entrainment of the oscillators through frequency and 
phase locking. However, in this model the coupling is assumed to be sufficiently weak that 
the amplitude is not affected. As a result, the model cannot describe amplitude instabili- 
ties. A system that does exhibit a rich variety of amplitude instabilities consists of coupled 
parametric oscillators and is the subject of this paper [0|§. 

There is a large literature on single parametric oscillators with periodically [[J or stochas- 
tically P modulated frequencies, perhaps also subject to thermal fluctuations J7|. A wide 
variety of linear and nonlinear, deterministic and stochastic, systems that exhibit energetic 
instabilities can be modeled as simple parametric oscillators undergoing the amplitude in- 
stability known as "parametric resonance" ||. It is therefore rather surprising that there 
is actually little work on coupled parametric oscillators. Work previous to our own that we 
are aware of is limited to a global parametric modulation that acts on all the oscillators 
in exactly the same way. Examples include parametrically pumped electrons in a Penning 
trap ||, pattern formation under global resonant forcing ||10|| , time-periodic loading of an 
elastic system , and nonlinear parametrically driven lattices [12| . 



In a previous paper we introduced a very simple, analytically solvable model of cou- 
pled parametric oscillators consisting of an infinite set of globally coupled harmonic oscil- 
lators subject to time-periodic piecewise-constant modulations with randomly distributed 
quenched phases ( "quenched" in this context means that the phase of the frequency modu- 
lation of each oscillator is set at time t = and then remains unchanged). We showed the 
appearance of a collective parametric instability: even though each individual oscillator is 
in its stable parameter domain, the average amplitude of the coupled system may diverge 
monotonically. This collective instability occurs when the amplitude of parametric mod- 
ulation is sufficiently large that the instantaneous frequency of the oscillators temporarily 
becomes imaginary. The instability is re-entrant with respect to the strength of the coupling 
of the oscillators and persists in the overdamped limit. In the presence of a saturating non- 
linearity M, it generates a pitchfork bifurcation, corresponding to a genuine second-order 
nonequilibrium phase transition (implying the spontaneous breaking of spatial symmetry 
and ergodicity). 

The purpose of this paper is to show that in addition to the instability described by 
monotonic growth of the mean amplitude, the globally coupled infinite system can also un- 
dergo transitions to a collective oscillatory instability with an intrinsic frequency that is not 
connected in a simple way to either the frequencies of the individual oscillators or to the 
frequency of the external modulation. A saturating nonlinearity turns this instability into 
a Hopf bifurcation generating a limit-cycle with the concomitant breaking of temporal sym- 
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metry and ergodicity [|]]. Furthermore, we show that even when the modulation amplitude 
is very small, monotonic growth can still occur but via a mechanism entirely different from 
that identified previously. 

In trying to explain these collective behaviors, one is impressed by the similarities with 
a simple model introduced earlier 0j of two coupled parametric oscillators modulated pe- 
riodically with a fixed phase difference 9 between their modulations. The behavior of this 
dimer when 9 = it is remarkably similar to that of the globally coupled infinite system, and 
the roots of the instabilities observed in the latter are already present in this simple system. 
At an even more primary level the seeds of these collective instabilities can be traced back 
to behaviors of single oscillators. In general an individual oscillator tends to synchronize 
with the external modulation whereas the coupling induces mutual synchronization between 
oscillators. These two synchronization tendencies cannot always be satisfied simultaneously. 
Coupling between oscillators can then be seen as leading to a sort of "selection" among the 
single oscillator modes, enhancing some (destabilization) and smoothing out others (stabi- 
lization) . 

In Section [TI| we introduce the globally coupled model. Section [TT1] is a review of the 
properties of a single linear parametric oscillator. Section [IV| establishes the mathematical 
setting for the analysis of the globally coupled system as a mean field problem, and typical 
numerical results are presented in Section |V[ In Section |VI| we collect these results in the form 
of a phase diagram that characterizes the behavior of the coupled system as the modulation 
parameters are varied. Here we discuss the boundaries between stable and unstable behavior 
and also between different instability regimes. Section [VI H offers a comparison between the 
mean field system and the dimer. Conclusions and a summary are presented in Section |V111. 



II. THE BASIC LINEAR MODEL 



Consider a set of TV unit mass linear parametric oscillators with displacements {xi}, each 
with a periodically modulated frequency and all harmonically coupled to one another. Here 
we restrict our analysis to coupled linear parametric oscillators. The nonlinear case will be 
presented elsewhere ||. The equation of motion of the i-th oscillator is given by 

k A 

Xi + 7^ + + ii{t)]xi = -— 22(xi - Xj) , (1) 

i=i 

with i = 1, • • • , N. 

Analytic results are possible with a simple piecewise-constant periodic modulation 
(period=T p ), defined as 

&(t) = A sgn[sinu p (t + n)] (2) 

where u p = 2n /T p , and the initial phase Tj is chosen at random for each oscillator. We are 
mainly interested in the mean amplitude 

1 N 

i=i 
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as a measure of the macroscopic behavior of the system. In the thermodynamic limit iV — > 
oo, this site average @ is equivalent to the average with respect to the random phase of 
the displacement of a single oscillator i, 

1 f Tp 

(x) = — / dr* , (4) 

which is independent of z. Equation ([]]) can then be reduced to a single mean-field differential 
equation, 

x + "fx + w%[l + £ T (t)]x = -k(x - {x)) (5) 

where we have dropped the index i and have now indexed the modulation with the initial 
phase chosen at random for each oscillator. Note that the average (x) must be evaluated 
self-consistently using Eqs. (§) and fl5|). 



III. REVIEW OF THE SINGLE PARAMETRIC OSCILLATOR 

The further analysis of the coupled system is clarified by first briefly examining the 
properties of a single parametric oscillator subject to the piecewise-constant frequency mod- 
ulation (fj) |p], |TT| , [T3^ [l~6[] , for which the equation of motion is simply 



x + i± + u [l + £ T (t)]x = 0. (6) 

The arbitrary phase r can be set to zero, but is retained explicitly in order to make clearer 
the connection with the coupled system. The equation of motion can be solved using Floquet 
theory (see the Appendix), and the result may be written as 

{t)+\x(t) ) = GAt) (i(0) + }x(0) ) (?) 



where G T (t) is the time-evolution operator, whose explicit expression is given in the Ap- 
pendix. 

Floquet theory provides information about the state of the system at the end of each 
modulation period, but it does not determine the time evolution within each period. If, as 
in most of the literature, stability conditions were our only concern, finding the Floquet 
eigenvalues would be sufficient. However, we are interested not only in the conditions giving 
rise to instabilities but also in the types of instabilities. Since exponential decay or growth 
of the oscillator amplitude is expected in any case, the time evolution may be conveniently 
investigated by Laplace transformation, 



x(s) = / e~ st x(t)dt 
Jo 

= [G T (s)] n x(0) + [G T (s)} 12 Qx(0) + x(0)) , (8) 
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where [G T (s)]jj is a matrix element of the Laplace transform of the time-evolution operator. 
In general, if x(s) has poles at Sj = Aj + iQj, the temporal behavior of the amplitude is 
expressed as 



x 



w=E 



Kit AQ.it 



cje J e 



(9) 



where the c/s are constants determined by the inverse Laplace transform and the initial 
conditions, 



Cj = 2ni lim (s — Sj / 

S— >Si 



G r (s)l x(0) + \G T (s)} (]-x{0)+x(0) 
J 11 L J 12 \2 



(10) 



The Laplace transform of the time-evolution operator is derived in the Appendix as 



GJs) 



l e -7T P /2 [ e (*+7)T P J_ G -l (Tp 



cosh 



R 



e- st G T (t)dt 



where I is the identity matrix and the resonance parameter R is defined as 



R = cos 



u + T p 



cos 



OJ-Tr. 



UJ 2 + + u 2 _ 



sm 



uj + T p 



sm 



with 



oj± — ujq \ 1 1 ± A 



7 

2uOr 



The poles of Eq. (|TT1) are determined by the condition 



cosh 



R = 



and their real parts are explicitly given by 



(11) 



(12) 



(13) 



(14) 



A, 



for 


R < -1 


for 


\R\ < 1 


for 


R > 1 , 



(15) 



while the imaginary parts are 



j - | for R < -1 



co p I j for R > 1 

Slf 1 

j ± — arccos(i?) for \R\ < 1 



UJ. n 



2tt 



(16) 



Here arccos(i?) lies in the range [0, it], j is an integer. Note that the real parts of the poles 
are in fact independent of j, and therefore one can drop the subscript j on A and rewrite 
the amplitude (|9]) as 



x(t) 



E 



Cj e 



(17) 
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Correspondingly, the relative weights of the different modes remain the same for all times: 
they all die (A < 0), or they all explode (A > 0), or they all keep their initial amplitudes 
(A = 0). None becomes relatively more dominant with time. 

Figure [l] illustrates the dependence of the poles on the resonance parameter. Consider 
first the real part A. When \R\ < 1 it is single-valued, R- independent, and negative (unless 
there is no damping 7, in which case it vanishes and the motion is purely oscillatory). At 
the bifurcation points R — ±1 the real parts become independent but remain negative until 
\R\ reaches the value 



which is greater than unity unless 7 = 0. Beyond \R\ = R c one of the A's becomes positive, 
leading to exponential growth of the oscillator amplitude. The condition \R\ = R c thus 
corresponds to the onset of parametric resonance or, more appropriately, of parametric 
instability. 

Consider next the imaginary parts Qj, shown in Fig. [I] for j = —2, ... 2 alternately as 
solid and dashed curves. The condition \R\ — 1 here also marks an interesting boundary at 
which the behavior of the imaginary parts changes from being /^-independent when \R\ > 1 
to independent when \R\ < 1. On the large- \R\ side of this boundary the frequencies Qj 
are simply proportional to the frequency u p of the modulation. When \R\ < 1, however, the 
oscillator frequencies change continuously with R and bear no simple relation to either u p 
or the natural frequency u . In Fig. ||we present the regions in parameter space (A, T p /T ) 
where \R\ > 1. The darker regions correspond to R > 1 and the lighter to R < —1. The 
stability boundaries \R\ — R c are also indicated by solid lines, the oscillator being unstable 
inside these boundaries. Note that the \R\ = R c and \R\ — 1 frontiers almost coincide 
because of the very low damping (7 = 0.1). As is well known, for small A the instability 
appears in the vicinity of T p « jT /2 (j integer). 

Although the boundary \R\ — R c is important in determining the transition from stable 
to unstable behavior for the single parametric oscillator, it does not play the same role for 
the coupled system, for which the boundary \R\ = 1 turns out to acquire further significance. 
Indeed, as we will see in the following section, each individual pole Sj = A+iQj gives rise to a 
collective pole in the coupled system, and these collective poles have different A's. Therefore, 
some modes may become unstable (A > 0), and therefore dominant, even when \R\ < 1, 
while others remain stable (A < 0). A very striking example is that of the single oscillator 
mode with a zero imaginary part, i. e., the one with j = when R > 1. If A is negative, it 
simply provides a monotonically decaying contribution to the oscillator displacement, and if 
A is positive it provides a monotonically growing contribution. However, its relative weight 
in the sum fllTD is extremely small and therefore this mode is practically not detectable in 
simulations of a single oscillator. Its contribution is overwhelmed by those of the oscillatory 
modes. However, in the coupled system this non-oscillatory mode may become dominant 
because it may become unstable while the oscillatory modes are still stable. This will lead 
to a monotonic explosion of the mean. 




(18) 
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IV. COLLECTIVE INSTABILITIES 



The mean field equation (0) can be rearranged as 

x + jx + {u$[l + Zr(t)] + k}x = k (x) (19) 
which describes a single parametric oscillator of frequency 



cu k = uj ox ll + A (20) 

driven by an "additive force" k (x) . Therefore, it can be solved with the time-evolution 
operator ( |A13|) of a single oscillator except that Eq. (|13|) must be replaced with the new 
shifted frequencies 



k ( 7 N 



\ 2uJ o 

The general solution of Eq. fll~9|) can be written as 



x(t) 
x(t) + lx{t) 



uj ± = uj QX jl±A+ — -[ ^- ) . (21) 



G ^ ( x(0)f: U) ) + kG ^ SI GAtTl ( <xS0> ) dt '- (22) 



Taking the average of Eq. (|22|) with respect to the random phase, we obtain a self-consistent 
equation for the mean amplitude, 

(x(t)) = ([G T (t)] n > x(0) + ([G T (t)] 12 ) (1 X (Q) + x(0)) + k f K(t - t) (x(t')) dt', (23) 

where the kernel is defined by 

K(t-t') = ([G T (t)G^(t')} 12 ) . (24) 



The kernel (J2J) depends only on the time difference because of the uniform distribution of 
the initial phases r. One can therefore solve the integral equation (|23|) by Laplace transfor- 
mation. The solution in Laplace space is given by 



[G r (s)]u)x(0) + ([GV(s)]i2> (fz(0) +x(0)) 

= ^ ~ • ( 25 ) 

l-A;([G T (s)] 12 



Equation (|25|) has a set of poles determined by the condition 



[G r (.s)] l2 ) = i, (26) 
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which differ from the ones derived for a single parametric oscillator. When k — > the right 
hand side of Eq. ( pq) diverges and the poles of the coupled system clearly approach those 
of the single oscillator. By explicitly performing the time integral in Eq. (|11|) and taking 
the average over the random phase r, one can obtain a rather complicated but explicit 
expression of Eq. (E6l) as 



cosh 



7T P 



-R 



x 



~2\ 2 

2 , T 



to 



~2\ 2 

7 



a; 



Axo> + (c<j 2 — u; 2 ) 2 



7 2 



sm 



oJ-T v 



cosh 



7 T p 



cos 



2 \2 



— kLU-(LU + — LU_) 



r 

4 



sm 



cosh I ^— - | — cos ' 



— ku)+U-y(uj + — uj_) sinh 



2 cosh 



7T P 



cos 



cos 



o;_T r . 



where 



(27) 



7 = 7 + 2s = 7 + 2A + 2iVL. 



(28) 



Finding the complex roots of Eq. ( |27| ) is difficult even numerically. Therefore, we first inves- 
tigate the collective modes graphically. Graphical inspection not only provides qualitative 
understanding of collective modes but also helps to identify suitable numerical algorithms. 

Coupling between the oscillators causes a frequency shift of each oscillator as indicated 
in Eq. fl2~iP . This is a trivial effect that is accounted for by calculating all single oscillator 
quantities such as R and the poles A + iQ using the shifted frequency. The interesting 
questions concern the way in which the couplings between the shifted oscillators affect the 
dynamics, and, in particular, whether instabilities and synchronization may be caused or 
suppressed by the coupling. A great deal can be learned about these dynamics starting 
from the single oscillator roots shown in Fig. [IJ. We thus investigate the coupled system in 
three typical regimes: (1) A single oscillator with the given (shifted) parameters is unstable 
(R < —R c ); (2) The single oscillator is stable and in the regime where its frequency is 
determined by the modulation frequency (1 < R < R c ); (3) The single oscillator is stable 
and its frequency bears no simple relation to the natural (shifted) oscillator frequency nor 
to the modulation frequency (\R\ < 1). The results are presented in the next section. 
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V. RESULTS 



To find the poles associated with collective modes we need to solve the coupled equations 

fcRe^[G T (s)]i 2 ^ - 1 = (29) 

Im([G T (s)] 12 ) = (30) 

for A and Q in s = A + iQ. In the next set of figures we present results for the three cases 
listed at the end of the last section. Each case is presented in two or three figures. 

In one set of figures we plot the left-hand sides of Eqs. ( P9"D and ( PU[ ) as contour lines 
in the space (A/u p , Q/lu p ) [the solid lines for Eq. (p9|), and the broken ones for Eq. (j30|)1. 
The parameters k = 1.28, ujq = 0.4, 7 = 0.01, and A = 1.0 are used in all the figures. 
The thick lines indicate zero contours. Therefore, the solutions we seek are the intersections 
of these two sets of thick lines, indicated by open circles. The poles of single oscillators 
with shifted frequency are indicated by black solid circles. There is of course an infinite 
number of solutions but we only exhibit those that fall within the scale of our figures. 
Associated with each case is also a figure representing a number of trajectories as a function 
of time that characterize the particular situation. For some of the examples we also show 
the associated phase space trajectories. Note that in this presentation we make a careful 
distinction between a single oscillator and an individual oscillator. The former refers to 
an independent oscillator with parameters and 7, while the latter refers to one of the 
oscillators in a coupled system with parameters Uk, 7, and k. 

Figure § shows the contours and poles for case (1). Here u>k = ui p /2, which corresponds 
to R = —1.0062. Since R < — 1, the poles of single oscillators appear on the lines Qj/u) p = 
{j — 2) each as a pair because in this regime there are two values of A associated with each 
Q. Only the poles for j = 1 are shown because the others (and there are an infinite number 
of them) are off the scale of the figure. Since for these parameters \R\ > R c = 1.00077, one 
of each pair of poles has a positive A. The amplitude of a single oscillator would therefore 
diverge exponentially. However, as the numerical calculations indicate, all the collective 
modes have A = —7/2 and thus the mean amplitude decays to zero despite the instability 
of single oscillators. 

Computer simulation results for various trajectories of a system of 100,000 oscillators 
with these parameters are shown in Fig. [|. We see that the mean amplitude is indeed zero, 

and that the deviation Ax = J (x 2 ) — (x) 2 diverges. The inset shows the same two trajec- 
tories as well as the diverging trajectory of an individual oscillator. Phase space trajectories 
of individual oscillators in the coupled system are shown in Fig. [5|. Each circle indicates a 
snapshot of an individual oscillator. Solid circles represent oscillators with positive modula- 
tion and open circles are ones with negative modulation at the time of the snapshot. Only 
1000 oscillators out of 100,000 are shown. 

The six snapshots show that with increasing time the phase volume increases (note 
the different scales in each snapshot), which is consistent with the divergent behavior of 
each oscillator and with the growth of the deviation Ax, and provide confirmation that 
there is indeed no mutual synchronization nor other kinds of organized collective motion. 
However, the persistent separation of solid and open circles into separate quadrants indicates 
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that individual oscillators are synchronized with the external modulation. Note that any 
individual oscillator moves clockwise, switching colors accordingly. However, since the phase 
of the modulation is chosen at random for each oscillator, at any time half of the oscillators 
have positive modulation and the other half negative, so that the pattern in preserved and 
an average over the phase gives (x) = 0. In this case, the effect of the coupling is only 
the frequency shift When the coupling strength is increased at fixed T p /Th, mutual 

synchronization will eventually be established. Then, the individual oscillators must lose 
synchronization with the modulation and consequently parametric instability is suppressed. 
Therefore, the coupling stabilizes the system in this parameter region. 

Interesting behavior is observed when u k = u p , that is, the single oscillator frequency is 
equal to the modulation frequency. This corresponds to case (2) on our list. The contours 
and poles for this case are shown in Fig. ^. Since R = 1.0001, the single oscillator has poles 
at Qj = ju p . Only the poles for j = are shown. This condition is close to parametric 
resonance, but since R is just below R c = 1.0031 all the single oscillator modes have negative 
A's and therefore single oscillator trajectories decay to the absorbing state x = 0. However, 
one of the collective modes has a pole with Q = and a positive A. All the other collective 
modes are dominated by this j = unstable mode and therefore the mean amplitude diverges 
monotonically. 

Computer simulation results for the associated trajectories are shown in Figure [?|. The 
mean decays slowly at the beginning and then diverges monotonically. The deviation di- 
verges as well, and does so more rapidly. The individual oscillator trajectory also diverges; in 
this example, although each single oscillator would be stable, the coupling causes individual 
oscillators in the system to become unstable. In other words, each oscillator is driven by the 
diverging mean on the right hand side of Eq. (0). The phase space trajectories in Fig. ||| 
show that after an initial transient (first three panels, where open circles are hiding most 
of the solid circles), individual oscillators in the coupled system oscillate with increasing 
amplitude about (x), while the mean (x) is moving away from the origin. In the long time 
limit, each oscillator "forgets" its initial conditions and it is driven by the mean. Therefore, 
while the phase of each individual oscillator is determined by the phase of the modulation, 
the amplitude in phase space is determined by the mean. Correspondingly, the oscillators 
become "amplitude-synchronized" through the mean. Until the synchronization is well es- 
tablished, the individual oscillators decay because the single oscillator modes have negative 
A. 

As in the previous case, the phases of all single-oscillator modes with Qj ^ become 
synchronized with the external modulation. The separation of solid and open circles in the 
three later-time panels in Fig. ^] indicates this synchronization. However, in contrast with 
the previous case, there is now a zero-frequency mode which does not have a phase to be 
synchronized. This mode is therefore not affected by the phase of each oscillator nor by 
the phase of the modulation. The zero-frequency mode shifts the center of oscillation away 
from x = in either the positive or the negative direction. In the presence of coupling the 
oscillators tend to follow the mean, and therefore shift in the same direction, thus breaking 
the system symmetry. Note that the open and solid circles no longer lie entirely in separate 
quadrants. It is precisely the excess of open circles relative to solid ones in the positive 
quadrants (and vice versa in the negative quadrants) that drives the motion of the mean 
further to the right. This mode can now dominate the behavior of the system because, 
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contrary to the situation of a single oscillator, each mode has a different A. In this example 
the effect of coupling is thus not only to shift the frequency according to (|2l|) , but also the 
more interesting collective symmetry-breaking monotonic divergence of the mean amplitude 
and the mutual synchronization of the individual oscillator amplitudes. 

Since the role of the coupling in this case is to suppress the single-oscillator modes with 
Qj 7^ and to allow the mode with Qj = to grow, this instability persists even for very 
large k. It should be noted that this collective instability is a different phenomenon from 
the monotonic growth reported previously 0. That instability occurs only for large A and 
disappears when the coupling strength k is increased (re-entrant transition). 

Finally, case (3) occurs when uj^ is between the two previous cases, so that \R\ < 1 and 
single oscillators are not in parametric resonance condition. However, Fig. || indicates that at 
least one of the collective modes has a positive A as well as non-zero Q. Therefore, the mean 
(x) oscillates with diverging amplitude. Recall that when \R\ < 1 the eigenfrequencies of 
single oscillators vary continuously with R but do not match the frequency of the modulation 
nor the natural frequency of the oscillator. Therefore, the individual synchronization to the 
external modulation plays no role and the phases of the oscillators are free to synchronize 
to one another through a synchronization to the phase of the mean. Computer simulations 
shown in Fig. ITU] confirm the oscillatory instability, and the phase space points of individual 
oscillators shown in Fig. [ll| testify to this phase synchronization. Although solid and open 
circles again form separate groups, the entire ring of open and solid circles alternates between 
the positive quadrants and the negative quadrants. Therefore, all oscillators are mutually 
synchronized. 

In the next section we collect the results for the coupled system into a phase diagram 
indicating regions of stability and instability of different types. 



VI. PHASE DIAGRAM 

A convenient way to summarize observations and characterize the instabilities system- 
atically is in the form of appropriate "phase diagrams" in which the stability boundaries 
are presented as a function of the system parameters. Since there are many parameters in 
this model the full diagram would involve a many-dimensional representation. Instead, we 
present the diagram in the two-dimensional space (A,T P /Tk) that characterizes the external 
modulation for a given set of oscillator parameters Uk, k, and 7. 

Figure O shows the phase diagram for the system parameters indicated in the caption. 



The colored regions denote unstable regimes, each color coding a particular type of instabil- 
ity. The yellow region denotes parameter ranges where the individual oscillators are unstable 
but the mean amplitude for the coupled system is zero ("incoherent unstable oscillations"). 
This region is associated with divergences of the first term in Eq. (^) but not of the second 
term. In this case the distinction between "single oscillators" and "individual oscillators" 
becomes moot, since the right hand side of Eq. (|19D plays no role. Oscillatory instabilities 
of the mean with positive A and nonzero Q are denoted in pink ("unstable spirals"). Blue 
regions ( " saddle nodes" ) denote monotonic divergence of the mean with one positive A and 
zero Q. Green regions ("unstable nodes") also denote monotonic instabilities but with two 
positive A's and zero Q. These latter three instabilities involve divergences of the second 
term in Eq. (E3). The first term may diverge as well, i.e., the A's of the single oscillator 
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and of the collective mode may be simultaneously positive. When this is the case we always 
observe the A of the collective mode to be larger than that of the single oscillator modes; 
the collective mode therefore dominates the dynamics. In establishing the terminology for 
various instabilities we have loosely followed the usual conventions of nonlinear dynamics. 

We can associate each specific case discussed in the previous section with a location 
on this sort of phase diagram. Thus the trajectories in Fig. |] are in a yellow regime of 
incoherent oscillatory divergence of each oscillator. Those of Fig. |7| are in a blue or green 
regime of monotonic divergence of the mean, and those of Fig. [10] are in the pink regime of 
unstable spirals. 

It is helpful to follow the behavior of the oscillator system across the various collective 
instability boundaries by considering the signs of Q and A as one increases the modulation 
amplitude A (thus moving upward vertically along the phase diagram) for different fixed 
values of the modulation period T p . Various associated bifurcation diagrams presenting A 
(solid lines) and Q (dotted lines) as a function of A are shown in Fig. [H| Consider first the 
period T p = 0.757)., shown in panel (a). As A increases, A changes sign, becoming positive 
at A c = 1.14 while Q remains positive. This represents a transition from a stable spiral to 
an unstable spiral (pink region in Fig. [12]). 

Consider next the period T p = 2.047)., shown in panel (b). Here A becomes positive at 
A c = 3.06 while Q remains positive. This therefore again marks a transition from a stable 
spiral to an unstable spiral (pink region). However, with a further increase in amplitude, Q 
eventually goes to zero at A' c = 3.51, where A bifurcates to two positive values (an unstable 
node, green region) via a saddle-node bifurcation. The oscillatory instability thus switches 
to a monotonic instability at this point. 

A different transition pattern is seen when T p = 3.07)., shown in panel (c). It begins 
with a stable spiral and switches to a stable node at A c = 2.91. With a further increase in 
amplitude the system undergoes a transition to a saddle node (blue region) at A' c = 3.04. 

A more complex transition pattern is shown in panel (d), in which the character of the 
instability changes several times along the T p = 3.97); line. As usual, at low amplitudes 
there is a stable spiral. At the point A c = 2.99 the system moves into an unstable spiral 
(pink region). The unstable spiral becomes an unstable node (very small green region in the 
phase diagram) via a saddle-node bifurcation at A' c = 3.30. A further transition to a saddle 
node (blue) occurs at A" = 3.4. 

We note that the phase diagram just described is quite rich and intricate. Instabilities 
cover even larger regions in parameter space, and do so with increasing intricacy, as the 
damping 7 decreases. A typical phase diagram for low damping is shown in Fig. [14. 

Increasing coupling or damping also leads to greater stability: the stability boundaries 
move "up" in the phase diagram when 7 and/or k are increased, indicating that a stronger 
modulation is needed to cause unstable behavior. Furthermore, the oscillatory instabilities 
eventually disappear with increasing modulation period T p , leaving only the monotonic 
collective instabilities. However, it should be noted that the monotonic instabilities for 
sufficiently large modulation period are simply due to an inversion of the effective harmonic 
potential and hence not due to any collective effects. That there is an inversion can already 
be suspected from the fact that at least one of the shifted frequencies in Eq. ( |2T[ ) becomes 
imaginary. 

An analysis of the system for large T p is fairly simple and instructive in elucidating the 
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source of instabilities more explicitly. In the adiabatic limit T p — > oo the single oscillator 
frequencies are frozen in time, half of them at frequency uj + and the other half at frequency 
£L, where 



The mean field equations of motion then are 

(£) + + 7 (x) + + ul (x) 



k 

%4 



± A. 



(31) 



{x 



(x}_ + 7 (x}_ + (x) 



2 x 7 ~ 



(32) 



where (■ ■ - )± indicates an average over the oscillators with the frequency u± respectively. 
This system can be diagonalized analytically in full generality. The eigenmodes of the 
coupled system are characterized by the complex frequencies 



nC 1 ) -1 _i_ 



k 



+ \ M 2 + 



k 
2^1 



n 



(2) 



1 + 



4^ 2 



- WA 2 



2ul 



(33) 



This case clearly illustrates the distinction between what we have called "single oscillator 
instabilities" and "chain instabilities". The former refer to the frequency ([31]) while the 
latter refer to (0). While the single oscillators would remain stable until A = 1 + k/2ujQ 
(at which point becomes imaginary), the chain becomes destabilized when A reaches 
the value \f\ + k/u^, where the imaginary part of becomes negative (note that the 
transition point is independent of 7). Beyond that the system is in a saddle-point /unstable- 
node instability region of non-oscillatory exponential growth due to potential inversion, 
which was the primary reason of instability in the previous work [0 . In the phase diagrams 
Figs. |T^| and this translates to a stability boundary that settles at A = 3 as T p — ► 00. 
Note that in particular the boundary A = yl + k/u^ remains valid in the overdamped limit 
7 — > 00, consistent with our early work |T7| . 

The origin of the instabilities presented as narrow blue and pink tongues in the low- A 
region of Fig. ^ is entirely different from the mechanism based on the temporarily inverted 
potential and unique to the underdamped case. In the following sections we will explain the 
cause of these instabilities using a dimer of parametric oscillators. 



VII. COMPARISON WITH PARAMETRIC OSCILLATOR DIMER 

In a recent paper || we reported results for a model of two coupled oscillators subject to 
parametric modulations with a phase difference 9. The equations of motion for this system 
are just the N = 2 version of Eq. ([!]): 
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k 

xi = -ujI[1 + £i(t)]xi - -(xi - x 2 ) - 7x1 

k 

x 2 = -o>o[l + 6(*)]a? 2 - 2 ( x 2 - »i) - 7x2- (34) 

The piecewise-constant periodic modulations of the two oscillators differ by a constant phase 
9, so that we can write Eq. (H) for this case as 

fi(t) = A sgn(sin(^ p t)) 

= A sgn (sin(^t + 9)) . (35) 

We want to investigate whether the mean position x = (xi + x 2 )/2 reproduces the macro- 
scopic behavior of the mean in the globally coupled model. 

In the absence of parametric modulations, a dimer has two eigenmodes: symmetric (or 
mutually synchronized, x 1 (t) = x 2 {t)) and antisymmetric (or mutually antisynchronized, 
Xi(t) = —x 2 (t)), with the former having a lower energy. When the parametric modulations 
are applied, these modes are in general no longer the eigenmodes of the dimer (except for 
9 = 0). However, the motions can be expressed as linear combinations of these modes and, in 
particular, the behavior of the mean of interest is reflected in the excitation of the symmetric 
mode by the parametric modulations. 

In previous work we focused on the behavior of the system as a function of the 
parameters A, T p = 2tc/u p , uq, 7, k, and 9. Among our conclusions is the fact that the regions 
of parametric instability are sensitively dependent on the phase difference 9. Of particular 
interest for the analysis in this paper is the behavior of the anti-phased dimer (9 = 71). 
This particular dimer captures many of the features of the mean field coupled system with 
unexpected accuracy. This assertion, which was originally based on our previous comparison 
of the regions of parametric resonance [Q], is reinforced when the dimer bifurcation diagrams 
are further refined to take into account details of qualitatively different trajectories, as we 
shall see below. 

For the piecewise constant parametric modulation (|35D the solution of the dimer problem 
is formally simple. The stability analysis is based on the eigenvalues {Xj} (|Ai| > |A 2 | > 
I A3 1 > IA4I) of the Floquet operator ||: parametric instability occurs when |Ai| > 1. To 
characterize different types of parametric instability we present bifurcation diagrams using 



the same color conventions as in the phase diagrams of the mean field model (Figs. [L2 



and IT3). If Im Ai 7^ 0, then clearly one has an oscillatory instability (pink regions). If 



Im Ai = 0, the instability can be either oscillatory or monotonic. Using the eigenvector 
corresponding to Ai as initial condition, we have determined whether or not x crosses zero 
during one period of the modulation. If it does, the point is assigned to a pink region. If 
it does not, the second largest eigenvalue A 2 will determine whether the point belongs to a 
blue (|A 2 | < 1) or green (|A 2 | > 1) region. The pink, blue and green regions are all caused 
by the instability of the symmetric mode. The yellow region, on the other hand, requires 
the instability of the antisymmetric mode and the decay of the symmetric mode. However, 
such a purely antisynchronous solution is forbidden in the 9 = it case @, which means that 
yellow regions cannot appear at all in the anti-phased dimer . 



Results for relatively large damping (7/cJo = 0.4) are presented in the bifurcation diagram 



of Fig. 15, which should be compared to Fig. O. The similarity between the two figures 
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is remarkable. Despite some extra green regions and the absence of the yellow tongue, 
one notices that the principal resonance regions of the dimer bifurcation diagram (1/2 < 
T p /Tk < 3/2) fits the same region in the mean field model almost exactly. The green regions 
connected to pink regions in the mean field model are also well mimicked by the dimer. 
Figure [16] illustrates the bifurcation diagram for a small value of damping (7/cJo = 0.01). 



Comparison with Fig. 14 shows that although the agreement between the two models is 



not as good as for higher values of 7/^0; the basic structure and similarities of the phase 
diagram and the bifurcation diagram are nonetheless preserved. The main differences are 
the complex pink patterns in region T v jT^ > 3 of the dimer. Also, as in the high-gamma 
case, the dimer has larger green regions than the mean field model, suggesting that the 
coupling in the latter plays a stronger role in stabilizing the system. But in spite of these 
differences and most importantly, the principal resonance region (1/2 < T p /Tk < 3/2) again 
shows an almost perfect match. 

In the dimer, competition between two kinds of synchronization plays a key role in the 
destabilization of the system: on the one hand, synchronization between each oscillator 
and its modulation; on the other hand, synchronization between the two oscillators. This 
competition is essentially governed by the values of A and k. Larger values of A favor the 
former, while larger values of k favor the latter. When the coupling is weak, the energy 
difference between symmetric and antisymmetric modes is small and both can be excited. 
In this case, the individual oscillators are nearly independent and the stability diagram 
of the dimer is similar to that of a single oscillator. As the coupling strength increases, 
the energy of antisymmetric oscillations increases until eventually only in-phase oscillations 
are energetically accessible. This mutually synchronized motion brings the system out of 
synchronization with the modulation. 

The stability diagram of the anti-phased dimer in the (Tp/T^, k) plane shown in Fig. |I7] 
illustrates the above story. First, consider T p /Tk « 0.5, where a single oscillator is in the 
main parametric instability region (R < —1, the first light grey region in the lower panel of 
Fig. |17|). For small k, the dimer is also unstable and still dominated by the antisymmetric 
mode (even though the symmetric mode cannot disappear, as mentioned above). However, 
as k increases, excitation of the antisymmetric mode becomes more difficult and the sym- 
metric mode becomes dominant, which stabilizes the system. When 0.5 < T p /Tk < la single 
oscillator is stable (\R\ < 1, between the dark and light grey regions). In this parameter 
region, individual oscillators do not have to be synchronous to the modulation (see Fig. p. 
They are free to become mutually synchronized and the system becomes unstable above 
a certain coupling strength. Since the symmetric mode dominates, this instability persists 
even for large k. Finally, when T p /T^ « 1, the single oscillator is again in an unstable region 
(R > 1, the dark grey region). Although the situation is similar to the first case, the individ- 
ual oscillators now have a zero-frequency mode which is not subject to synchronization with 
the modulation. Therefore, the zero-frequency mode of the two oscillators can be mutually 
synchronized, which produces a monotonic growth of the mean. Again, this instability does 
not have strong dependence on the coupling strength and persists even for large k. 

The good agreement between the dimer and mean field models is not merely a coinci- 
dence. Consider one particular test oscillator in the globally coupled system. Now define 
the group of its "friends" as the set of all oscillators whose modulation phase lies within 
an interval of ±7r/2 around its phase. The group of its "enemies" comprises the set of all 
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other oscillators, whose average modulation phase is opposite to that of the test oscillator. 
When the oscillators are synchronized to the modulations, "friends" are also mutually syn- 
chronized to one another, regardless of the coupling. When the coupling increases, there is a 
competition between two kinds of synchronization: synchronization between the "enemies" 
and "friends" , and synchronization with their own modulations. This situation is similar to 
the anti-phased dimer, which helps explain the remarkable similarity between the instability 
diagrams of the dimer and mean field models. Notice that this description becomes exact in 
the quenched limit, where the "friends" and "enemies" can be represented respectively by 
(x), and (x)_. From this perspective, the choice 9 = n for the dimer appears as a natural 
one, being not only a particularly symmetric case in the general dimer problem H but also 
the effective phase difference between the two groups. 

VIII. DISCUSSION AND CONCLUSIONS 

We have investigated collective instabilities of an infinite set of globally coupled linear 
oscillators driven by time-periodic piecewise linear modulations with random initial phases. 
These instabilities occur in certain parameter regimes (and not in others), and we have 
produced phase diagrams as a function of system parameters indicating detailed stability 
boundaries and types of instabilities. 

Instabilities arise from phase synchronization (although not all phase synchronization 
leads to instability), and there are two possible competing synchronization mechanisms: 
synchronization of individual oscillators with the external modulation ("modulation syn- 
chronization"), and mutual synchronization between oscillators ("mutual synchronization"). 
In the absence of the external modulation, only mutual synchronization is possible. On the 
other hand, when the coupling is absent, a single oscillator in parametric resonance condition 
(\R\ > 1) synchronizes with the external modulation. Note that parametric instability of a 
single oscillator requires modulation synchronization, and in general any instability requires 
an external modulation to pump energy into the system. 

In the presence of external modulation with random phases and coupling, even in a pa- 
rameter regime where either alone would lead to synchronization, it is not possible for both 
types of synchronization to occur simultaneously. Intuition might lead to the conclusion that 
sufficiently strong coupling would favor mutual synchronization thereby destroying modu- 
lation synchronization (and hence stabilizing the system). Intuition might also lead one to 
conclude that weak coupling necessarily results in a dominance of modulation synchroniza- 
tion and hence to instabilities of individual oscillators not related to one another. However, 
we have shown that this intuitive picture would be incomplete because it does not account 
for the existence of collective parametric instabilities due to the coupling and entirely absent 
in a single parametric oscillator. Here we discuss our results in terms of the two competing 
synchronization mechanisms. 

In a globally coupled model, the individual oscillators are coupled to the mean (x), and 
thus mutual synchronization can be thought of as synchronization between the individuals 
and the mean. Note that there is no mutual synchronization when (x) = (because if there 
were, the mean would not be zero). 

Now suppose we are in the parameter regime R < — 1 where individual uncoupled oscil- 
lators are parametrically unstable. When these oscillators are coupled, one can imagine one 
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of two possible scenarios. If the coupling leads to mutual synchronization, the individual 
oscillators can no longer be synchronous with the external modulation and therefore the 
coupled system has been stabilized by the coupling. On the other hand, if the coupling does 
not lead to mutual synchronization but instead there is modulation synchronization, then 
the oscillators may be individually unstable but with (x) = 0. Our results show that the 
second scenario is the correct one for sufficiently small values of k, as shown in the yellow 
regions of incoherent instability in Fig. [TJ|. Modulation synchronization has "won" . On the 
other hand, there is a coupling energy cost to the lack of mutual synchronization, which 
slows down the instability of individual oscillators relative to their uncoupled amplitude 
growth. For larger k the first scenario takes over, and the yellow region disappears above a 
certain value of the coupling. 

Next suppose that we are in the other parametric resonance regime, R > 1. The situation 
is in some ways similar to the previous case but there is a major difference: there is now 
a mode, the j = mode, whose frequency is zero and which therefore (contrary to the 
other modes) need not (indeed cannot) synchronize to the modulation. The amplitude of 
this mode can grow monotonically in either direction, and the coupling among oscillators 
leads to a tendency for the zero-frequency mode of all the individual oscillators to move in 
the same direction. Thus while the growth rate of the j ^ modes is reduced by coupling 
due to the lack of mutual synchronization, that of the j = mode is enhanced because the 
coupling fosters mutual synchronization of this mode. 

If the individual oscillators are not in regimes of parametric instability (\R\ < 1), there 
is no synchronization to the modulation and the oscillators are free to synchronize with 
one another. Mutual synchronization is thus established and the mean (x) becomes oscilla- 
tory with the same frequency as that of individual oscillators in the coupled system. The 
oscillatory mean drives the system into unstable states via the mean field coupling. 

We have found that the instability bifurcation diagram for a simple anti-phased dimer 
model reproduces the phase diagram for the mean field system with surprising accuracy. 

In the absence of external modulation, a dimer has two eigenmodes: symmetric or mutu- 
ally synchronized (lower energy), and antisymmetric or mutually anti-synchronized (higher 
energy). In the presence of time-periodic piecewise linear modulations which are exactly out 
of phase on the two coupled oscillators (9 — 7r), a competition between these two modes 
(which are no longer eigenmodes) ensues. This competition is in many ways similar to the 
competition between modulation synchronization and mutual synchronization described for 
the globally coupled system, and here again it determines the stability of a dimer. When 
the coupling is weak, the energy difference between symmetric and antisymmetric modes is 
small and both can be excited. In this case, the individual oscillators are nearly indepen- 
dent and the stability diagram of the dimer is similar to that of a single oscillator. The 
synchronization of each oscillator with its modulation dominates the behavior, and insta- 
bilities thus represent boundless excitation of the antisymmetric mode. With increasing 
coupling the energy of the antisymmetric mode increases until it is too high to be excited. 
Only the symmetric mode can be excited, i.e., the oscillators become mutually synchronized. 
The synchronization with the modulation is thus destroyed and the associated parametric 
instability is suppressed. 

Although the similarity between the dimer and globally coupled models is remarkable, 
they also exhibit various important differences. In the dimer model, mutual synchronization 
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involves only two oscillators. On the other hand, in the global coupling model an oscillator 
must be synchronous with essentially all others to create collective motion. Therefore, in 
the thermodynamic limit iV — > oo, the collective instability in the globally coupled system 
is a genuine phase transition, whereas the instabilities in the dimer are simple bifurcations. 
Nevertheless, the stability boundaries and dynamics of the mean amplitudes in both cases 
show impressive similarities. 

An interesting case that is in some sense "in between" these two and that promises in- 
teresting new features is that of a one-dimensional chain of oscillators with nearest-neighbor 
coupling. When the phase of the modulations of the oscillators in the chain is chosen at 
random, there is a significant chance that both neighbors of any given oscillator have a 
phase "similar" (suitably defined within some range) to its phase. In this case, the middle 
oscillator can easily establish simultaneously both mutual synchronization with its neighbors 
and synchronization with the modulation. Therefore, locally this oscillator could become 
unstable. On the other hand, if the neighbors of a given oscillator are modulated with phases 
opposite to its own modulation phase, that oscillator may be stabilized. Therefore, the spa- 
tial pattern of the modulation phase is expected to play an important role and instability 
may become wavelength dependent, suggesting spatial pattern formation. Such patterns can 
of course not be observed in either a dimer or a globally coupled model. A detailed analysis 
of such systems will be presented elsewhere j| . 
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APPENDIX A: TIME-EVOLUTION OPERATOR OF A SINGLE PARAMETRIC 

OSCILLATOR 

We solve the equation of motion (|j) using a standard Floquet method. The damping 
term can be eliminated by introducing a new variable y defined as 

x(t) = e-^ 2 y(t). (Al) 

Equation (|5]) becomes 

y + u 2 {t)y = (A2) 

with the time-dependent frequency 

u\t)^ul\l + Ut)}-l 2 /^ (A3) 

The solution of the undamped frequency-modulated oscillator (|A2|) can be expressed in 
terms of the time-evolution operator from t = to t, g T (t), as 

^-^(ym (A4) 
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For a piecewise constant modulation such as Eq. (0), the explicit form of the time-evolution 
operator is known. Using its periodicity and composition property we note that if t = nT + u 
then 



g T (t) = g T (nT + u) = g T (u)g T (nT) = g T (u) [g T (T)} n . 

It is thus sufficient to find g T {t) for t E [0, T] 

When the phase is r E [0, T/2], the frequency varies as 



(A5) 



while for r E [T/2,T], 
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where the u)± are defined in Eq. (|T|). During each constant frequency time window the 
system evolves according to the propagator of a simple harmonic oscillator of the appropriate 
frequency. This propagator is well known: 



g± (*,*') 



cos[w±(t-f)] at sin[o; ± 
-u± sin[a>±(t — t')] cos[uj±(t — £')] 
The full operator g T (t) can be expressed as products of the g±. For the case 
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^g + (t,T-r)g_(T-r,^-r)g + (^-r,0) for t G [T — r, T] 



and for the case (|A7|) , 
'g-M) 
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These expressions can be further simplifed by using time-translation symmetry, g±(t + u, t' + 
u) = g±(t,t') = g±(t — t',0). We can thus simplify our notation and set g±(t, 0) = g±(t) 
with no loss of information. 
In particular, for t = T 

g r (T) = S T g (T)S; 1 (All) 

where 

goCO = g_ (£) g+ K) , S T ^g + (r), for rG[0,|] 

V J V J (A12) 
g (T) = g + (|) g_ (|) , S r = g_(r-^), for rG[|,T]. 

Since only S r depends on r, the trace, determinant, and eigenvalues of g r (T) do not depend 
on the random phase. 

Transforming back to the original variables, we finally obtain the time-evolution operator 
for x(t) and x(t), 

G T (t) = e-^ 2 g T (t). (A13) 

Next we derive an explicit expresion for the Laplace transform of the evolution operator 
with the help of Eq. 
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(A14) 



G T (s) = / e- st G T (t)dt 
Jo 

= ^[e- sT G T (T)] m / e- st G T {t)dt 

rn=0 ^° 

= (l-e- sT G T (T)y 1 [ e~ st G T (t)dt 

Jo 

i e - 7 T/ 2 [ e ( s +7)Tj _ G r (T)] e~ st G T (t)dt 
cosh[(s + 5) T] -i? 
where 7 is an identity matrix and 
i? = I^ T / 2 Tr G (T) 

= ^Tr go(T) 

= cos(-^-) cos(^-) - 2 sin(— ) sin(— ) (A15) 

The geometric series in Eq. ( |A14 ) converges only when 

|| e - sT G (T)||<l. (A16) 
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FIGURES 
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FIG. 1. Real (upper graph) and imaginary (lower graph) parts of the poles Sj as a function of 
R with j/uj p = 0.1. The behavior associated with these poles is described in detail in the text. 
Note the significance of the boundaries R = ±1, and of the boundaries R = ±R C where one A 
becomes positive. The latter mark the onset of single-oscillator parametric instability. 




boundaries \R\ = R c . The damping 7 = 0.01 and the frequency loq = 0.4. 
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A/ Wp 

FIG. 3. Contour plots of the left-hand side of Eq. ( ]29| ) (solid lines) and Eq. ( pop (broken lines) 
for uj = 0.4, k = 1.28, 7 = 0.01, A = 1.0, and T k = 2T p , which lead to R = 1.00617. Poles of 
a single oscillator are indicated by solid circles. Thick lines correspond to solutions of Eqs. ( |29| ) 
(solid) and ([}(]) (broken). The intersections of the solid and dotted thick lines, indicated by open 
circles, are thus solutions of Eq. (p6|). Although there is a single oscillator pole with positive A 
there is no positive-A collective mode. 




23 




24 



-0.01 0,01 

FIG. 6. Contour plots of the left-hand side of Eqs. ( p9| ) (solid lines) and Eq. ( p0| ) (broken lines) 
for to = 0.4, k = 1.28, 7 = 0.01, A = 1.0, and T k = T p , which lead to R = 1.0001. Poles of 
a single oscillator are indicated by solid circles. Thick lines correspond to solutions of Eqs. ( |29| ) 
(solid) and ([}(]) (broken). The intersections of the solid and dotted thick lines, indicated by open 
circles, are thus solutions of Eq. (|26|). Although the single oscillator poles have all negative A, 
the collective modes include a pole with positive A and S7 = 0. This mode diverges exponentially 
without oscillation. 
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FIG. 7. Trajectories associated with Fig. |6|. Thick solid line: mean (x). Thin solid line: 
deviation Ax. The inset also includes the trajectory x of an individual oscillator (broken line). 
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FIG. 9. Contour plots of the left-hand side of Eqs. (|29j ) (solid lines) and Eq. ( |30|) (broken lines) 
for u = 0.4, k = 1.28, 7 = 0.01, A = 1.0, and T k = 4T p /3, which lead to R = -0.01037. Poles 



of a single oscillator are indicated by solid circles. Thick lines correspond to solutions of Eqs. (29) 
(solid) and (^) (broken). The intersections of the solid and dotted thick lines, indicated by open 
circles, are thus solutions of Eq. (^(j). Although the single oscillator poles all have negative A, the 
collective modes include a pole with positive A and nonzero. This mode diverges exponentially 
with oscillation. 
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FIG. 10. Trajectories associated with Fig. || Thick solid line: mean (x). Thin solid line: 
deviation Ax. The inset also includes the trajectory x of an individual oscillator (broken line). 
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FIG. 11. Phase point snapshots of 2000 individual oscillators in the system associated with 
Figs. P and 10. Note the scale changes with increasing time. 
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FIG. 12. Phase diagram for the mean field model with oscillator parameters ojq = 0.4, k = 1.28 
and 7 = 0.16. White regions denote stable regimes. The various instability regimes are color coded 
as indicated. The characteristic behavior in each instability regime is described in the text. 
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FIG. 13. Bifurcation diagrams showing A (solid lines) and 17 (dotted lines) with changing 
modulation amplitude for various values of the modulation period. Panel (a): T p /T k = 0.75; (b): 
T p /T k = 2.04; (c): T p /T k = 3.0; (d): T p /T k = 3.9. The behavior implied by these diagrams is 
discussed in detail in the text. 
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FIG. 14. Phase diagram for the mean field model with oscillator parameters ojq = 0.4, k = 1.28 
and 7 = 0.01. 
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FIG. 15. Bifurcation diagram for the anti-phased dimer with the same parameters as in Fig. ^ 
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FIG. 17. Upper panel: Bifurcation diagram for the anti-phased dimer in the (T p /Tk,k) plane 
for A = 0.9, 7 = 0.01, and ojq = 0.4, using the same color convention as in Fig. Lower panel: 
R (Eq. 013)) as a function of Tp/T^; light grey areas denote R < — 1 and dark grey areas -R > 1. 
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